############################################################################################
# Replication code and data for the manuscript:
# How do coalition governments legislate?
# Drafting and negotiating bills in the executive phase of law-making 
############################################################################################
# Date 30.7.2020


# General comment:
############################
#  We have included only final variables employed in the models, not the raw data.
#  However, the raw dataset is also available, it contains all the variables we were able to collect.
#  Variables are typically transformed using natural log. 
#  For raw data, see the unprocessed dataset. We have modified variables by a separate script (not shown here).
#  This script has been modified during the review process. Some variables used in previous models were relaxed from the
#  final model of the final manuscript.

# Script is divided into several sections:
#  1) Install packages
#  2) Import dataset
#  3) Subsetting of variables
#  4) Standardization of variables
#  5) Descriptive statistics
#  6) Main regression model - poisson regression
#  7) regression diagnostics
#  8) plots of predicted values ŷ based on values of x
#  9) means plots for categorical variables

# Required packages - must be installed
library(lme4)
library(texreg)
library(haven)
library(ggplot2)
library(ggeffects)
library(tidyverse)
library(corrplot)
library(car)
library(tidyverse)
library(hrbrthemes)
library(viridis)
library(ggpubr)
library(AER)
library(car)
library(boot)
library(pastecs)
library(stargazer)

# Setting working directory to your data file location
setwd("to your source file location")   #copy-paste the address of the folder

setwd("~/") #or to the source file location

# Importing dataset - SPSS file
BILL_CHANGE_DATASET <- read_sav("BILL_CHANGE_DATASET.sav")

#as STATA file (optional)
# BILL_CHANGE_DATASET <- read_dta("BILL_CHANGE_DATASET.dta")


# Hypothesis and variables: 
############################


# H1: The higher the policy distance on a bill to the political party holding the ministry 
# submitting the bill to the coalition compromise, the more changes the bill undergoes in the executive phase.

# H2: The higher the saliency of a bill to other coalition parties (than the party holding the ministry 
# submitting the bill), the more changes the bill undergoes in the executive phase.  

# H3a: If the government bill is listed in the document that presents the priorities of the
# coalition government, the bill undergoes less changes in the executive phase.  

# H3b: The more the initiating ministry changes the government bill before it officially 
# enters the executive phase, the less changes the bill undergoes in the executive phase.

# H4: The more closely the government bill is scrutinized by institutions reviewing its 
# regulatory and legal quality, the more changes the bill undergoes in the executive phase.  




# Control hypotheses:
#	C: A bill that amends the existing laws (amendment) is changed less than a self-standing bill (new law). 
#	C: A bill initiated closer to the end of the election period is changed less.  
#	C: Bills that stay in the executive phase longer probably face difficulties and tend to have more alterations. 
#	C: A longer (number of words) bill is changed less.  

#	Name of the variable in the dataset	-	Label

#	H1 Distance_COL_COMP_s	-	Distance from coalition compromise - Martin and Vanberg (2014) based on ChES DATASET
#	H2 SALIENCY_COAL_im_mean - "Saliency of a bill to other coalition parties (computed as the mean)
#	H3a -  Legislative_plan	-	Government's Legislative Plan - dummy variable
# H3b - "Draft_DIFF_BI" - Ministerial phase change (control)
#	H4 - "LRV_meeting"	-	Number of LCG meetings


#	Control - "Amended_law"	-	Amendment  - dummy variable
#	Control - "Draft_lenght_ln"	-	Bill´s length (N words logged) - number of words divided by 1000 and logged
#	Control - "Days_gov_approval_ln"	-	Days between the final gov. approval and the end of IMC (ministerial phase) (N logged) 
#	Control - "Days_to_elections_ln"	-	Days to next parliamentary elections (N days logged)


# Control dummies
# Government - 2 Sobotka government, 4 Nečas governments (as.factor command produce dummy)
# Proposer_cat - string variable indicating the ministerial portfolio
# Dimenze_Chapel_Hill_OTHER_CAT - CHES dataset, including imputed cathegory "other".
# PolicyArea - detailed policy area (own coding based on the description of the bill)

# Excluded variables (and relaxed hypotheses) after revision of the manuscript:

#	H1: A bill initiated by the ministry held by the main coalition party is changed more. 
#	H2: The higher the number of disagreements in the comment procedure relating to a bill
#     from ministries held by other coalition parties, the more changes the bill undergoes. 
#	H3: The higher the saliency of a bill to the political party holding the ministry submitting 
#     the bill, the less changes the bill undergoes 
# H5: The larger the policy distance on a bill to the political party holding the ministry
#     submitting the bill to the median position of the legislative party, the more changes the bill undergoes.


#	H1 Main_ruling_party_DB	-	Drafter: Main ruling party - dummy variable
#	H2 Comments_coal_disputed_ln	-	Disagreements comments - coalition partner (N logged)  
#	H3 Saliency_CH_H	-	Saliency of a bill to the party - Zubek and Klüver (2015). 

# H5 Distance_LEG_MEDIAN -	Distance from legislative median party position - Martin and Vanberg (2014) based on ChES DATASET
#	Control - "Complexity"	-	Complexity (N of policy areas) - number of policy areas stated in EKLEP
#	Control - "EU_transposition"	-	EU transposition - dummy variable if stated in EKLEP



#Subseting datasets - run before analysis!
###########################################################################################
#Datatset DRAFT_BILLS = the unit of analysis is a bill draft

# Selecting variables for the Table 1
TABLE_1 <- BILL_CHANGE_DATASET[, c("DRAFT_DIFF_PAR_SENAT_BI", "DRAFT_DIFF_IMC_PRINT_BI")]


# Selecting only variables for the regression models and clearing missing data!
BILL_CHANGE_DATASET <- BILL_CHANGE_DATASET[, c("ID", "Draft_DIFF_UNI",
                                               "Draft_DIFF_BI",
                                               "DRAFT_DIFF_GOV_PAR_UNI",
                                               "DRAFT_DIFF_GOV_PAR_BI",
                                               "DRAFT_DIFF_IMC_PAR_UNI",
                                               "DRAFT_DIFF_IMC_PAR_BI",
                                               "Draft_DIFF", "DRAFT_DIFF_GOV_PAR","DRAFT_DIFF_IMC_PAR", "Complexity",
                                             "Main_ruling_party_DB", "Comments_coal_partner_ln", "Comments_same_party_ln", 
                                             "Comments_same_party_disputed_ln",
                                             "Comments_coal_disputed_ln","COM_IMP_DB",
                                             "Saliency_CH_H", "Distance_COL_COMP_s",
                                             "Distance_LEG_MEAN", "Distance_LEG_MEDIAN",
                                             "LRV_meeting", "Legislative_plan", 
                                             "Days_to_elections_ln","Amended_law", "EU_transposition",
                                             "Draft_lenght_ln", "Days_gov_approval_ln",
                                             "Proposer_cat", "Government", 
                                             "Dimenze_Chapel_Hill_OTHER_CAT",
                                             "PolicyArea", 
                                             "Comments_coal_partner_disputed", "Days_to_elections",
                                             "Days_gov_approval", "Draft_lenght", "SALIENCY_COAL_im_zero",
                                             "SALIENCY_COAL_im_mean")]
BILL_CHANGE_DATASET <-  na.omit(BILL_CHANGE_DATASET)


#Unstandardized dataset
BILL_CHANGE_DATASET_UNSTANDARDIZED  <- BILL_CHANGE_DATASET

# Selecting government - for specific analysis not presented in manuscript
NECAS <- BILL_CHANGE_DATASET[ which(BILL_CHANGE_DATASET$Government == 2), ]
SOBOTKA <- BILL_CHANGE_DATASET[ which(BILL_CHANGE_DATASET$Government == 4), ]


#STANDARDISATION - Gelman method 2sd
##############################################################################
# do not run this part of the script if you want unstandardized results

twoSD <- function (x) {
  (x - mean(x, na.rm = TRUE))/(2 * sd(x, na.rm = TRUE))
}

BILL_CHANGE_DATASET$Amended_law <- twoSD(BILL_CHANGE_DATASET$Amended_law)
BILL_CHANGE_DATASET$EU_transposition <- twoSD(BILL_CHANGE_DATASET$EU_transposition)
BILL_CHANGE_DATASET$Draft_lenght_ln <- twoSD(BILL_CHANGE_DATASET$Draft_lenght_ln)
BILL_CHANGE_DATASET$Legislative_plan <- twoSD(BILL_CHANGE_DATASET$Legislative_plan)
BILL_CHANGE_DATASET$Main_ruling_party_DB <- twoSD(BILL_CHANGE_DATASET$Main_ruling_party_DB)
BILL_CHANGE_DATASET$Comments_coal_partner_ln <- twoSD(BILL_CHANGE_DATASET$Comments_coal_partner_ln)
BILL_CHANGE_DATASET$Comments_same_party_ln <- twoSD(BILL_CHANGE_DATASET$Comments_same_party_ln)
BILL_CHANGE_DATASET$Comments_same_party_disputed_ln <- twoSD(BILL_CHANGE_DATASET$Comments_same_party_disputed_ln)
BILL_CHANGE_DATASET$Comments_coal_disputed_ln <- twoSD(BILL_CHANGE_DATASET$Comments_coal_disputed_ln)
BILL_CHANGE_DATASET$Days_gov_approval_ln <- twoSD(BILL_CHANGE_DATASET$Days_gov_approval_ln)
BILL_CHANGE_DATASET$COM_IMP_DB <- twoSD(BILL_CHANGE_DATASET$COM_IMP_DB)
BILL_CHANGE_DATASET$Saliency_CH_H <- twoSD(BILL_CHANGE_DATASET$Saliency_CH_H)
BILL_CHANGE_DATASET$Distance_COL_COMP_s <- twoSD(BILL_CHANGE_DATASET$Distance_COL_COMP_s)
BILL_CHANGE_DATASET$Distance_LEG_MEAN <- twoSD(BILL_CHANGE_DATASET$Distance_LEG_MEAN)
BILL_CHANGE_DATASET$Distance_LEG_MEDIAN <- twoSD(BILL_CHANGE_DATASET$Distance_LEG_MEDIAN)
BILL_CHANGE_DATASET$Complexity <- twoSD(BILL_CHANGE_DATASET$Complexity)
BILL_CHANGE_DATASET$Days_to_elections_ln <- twoSD(BILL_CHANGE_DATASET$Days_to_elections_ln)
BILL_CHANGE_DATASET$LRV_meeting <- twoSD(BILL_CHANGE_DATASET$LRV_meeting)
BILL_CHANGE_DATASET$SALIENCY_COAL_im_mean <- twoSD(BILL_CHANGE_DATASET$SALIENCY_COAL_im_mean)



BILL_CHANGE_DATASET_MINISTERIAL <- BILL_CHANGE_DATASET

#if as independent variable in the model

BILL_CHANGE_DATASET$Draft_DIFF_BI_UNSTANDARDIZED <- BILL_CHANGE_DATASET$Draft_DIFF_BI

BILL_CHANGE_DATASET$Draft_DIFF_BI <- twoSD(BILL_CHANGE_DATASET$Draft_DIFF_BI)




# Creating dummy variables (as.factor does not function with all the packages e.g. interactions)
BILL_CHANGE_DATASET$SOBOTKA_GOV <- ifelse(BILL_CHANGE_DATASET$Government == '4', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_1 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '1', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_2 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '2', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_3 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '3', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_4 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '4', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_5 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '5', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_6 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '6', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_7 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '7', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_8 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '8', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_9 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '9', 1, 0)
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_10 <- ifelse(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT == '10', 1, 0)

# saving as dta file for STATA software - optional
write_dta(BILL_CHANGE_DATASET, "BILL_CHANGE_STANDARDIZED.dta")


# Descriptive statistics
#############################################################################################################

# Histograms of the DV

# Ministerial phase
ggplot(data = BILL_CHANGE_DATASET_UNSTANDARDIZED) + 
  geom_histogram(mapping = aes(x = Draft_DIFF_BI), binwidth = 5) +               
  labs(x = "Ministerial phase change", y = "Frequency")                       

# Cabinet phase - DV
ggplot(data = BILL_CHANGE_DATASET) + 
  geom_histogram(mapping = aes(x = DRAFT_DIFF_GOV_PAR_BI), binwidth = 5) +               
  labs(x = "Cabinet phase change", y = "Frequency")                       

# Total change in the executive phase
ggplot(data = BILL_CHANGE_DATASET) + 
  geom_histogram(mapping = aes(x = DRAFT_DIFF_IMC_PAR_BI), binwidth = 5) +               
  labs(x = "Executive phase change", y = "Frequency")                       


# Descriptive statistics 

# Table 1

library(stargazer)
stargazer::stargazer(BILL_CHANGE_DATASET_UNSTANDARDIZED[c("Draft_DIFF_BI","DRAFT_DIFF_GOV_PAR_BI","DRAFT_DIFF_IMC_PAR_BI")], type = "text",
          title="Table 1: Bills change in different phases of the legislative process (in percent)", digits=1, out="table1.txt",
          covariate.labels=c("Ministerial phase","DV: Cabinet phase","Executive phase (total)"))

stargazer::stargazer(TABLE_1, type = "text", summary=TRUE, min.max=TRUE, mean.sd = TRUE, 
          nobs = FALSE, median = FALSE, iqr = FALSE, title="Descriptive statistics", digits=1, out="table1.txt")

stargazer::stargazer(TABLE_1, type = "text")

library(pastecs)
options(scipen=100)
options(digits=3)
scores <- cbind(BILL_CHANGE_DATASET$Draft_DIFF_BI,
                BILL_CHANGE_DATASET$DRAFT_DIFF_GOV_PAR_BI,
                BILL_CHANGE_DATASET$DRAFT_DIFF_IMC_PAR_BI)
stat.desc(scores)

scores <- cbind(TABLE_1$DRAFT_DIFF_PAR_SENAT_BI,
                TABLE_1$DRAFT_DIFF_IMC_PRINT_BI)
stat.desc(scores)

# Table 2

stargazer(mydata[c("mpg","hp","drat")], type = "text",
          title="Descriptive statistics/selected variables", digits=1, out="table2.txt", flip=TRUE,
          covariate.labels=c("Miles/(US)gallon","Gross horsepower","Rear axle ratio"))

# Unstandardized values
scores <- cbind(BILL_CHANGE_DATASET_UNSTANDARDIZED$Distance_COL_COMP_s,
                BILL_CHANGE_DATASET_UNSTANDARDIZED$SALIENCY_COAL_im_mean,
                BILL_CHANGE_DATASET_UNSTANDARDIZED$Legislative_plan,
                ILL_CHANGE_DATASET_UNSTANDARDIZED$LRV_meeting,
                BILL_CHANGE_DATASET_UNSTANDARDIZED$Amended_law,
                BILL_CHANGE_DATASET_UNSTANDARDIZED$Days_to_elections,
                BILL_CHANGE_DATASET_UNSTANDARDIZED$Days_gov_approval,
                BILL_CHANGE_DATASET_UNSTANDARDIZED$Draft_lenght)
stat.desc(scores)



# Standardized values - used in the model
scores <- cbind(BILL_CHANGE_DATASET$Distance_COL_COMP_s,
                BILL_CHANGE_DATASET$SALIENCY_COAL_im_mean,
                BILL_CHANGE_DATASET$Legislative_plan,
                BILL_CHANGE_DATASET$LRV_meeting,
                BILL_CHANGE_DATASET$Amended_law,
                BILL_CHANGE_DATASET$Days_to_elections,
                BILL_CHANGE_DATASET$Days_gov_approval,
                BILL_CHANGE_DATASET$Draft_lenght)
stat.desc(scores)


# Correlation matrix
###################################
COR_MATRIX <- BILL_CHANGE_DATASET[, c("DRAFT_DIFF_GOV_PAR_BI", "Distance_COL_COMP_s",
                                      "SALIENCY_COAL_im_mean","Legislative_plan", "Draft_DIFF_BI",
                                      "LRV_meeting","Amended_law", "Days_to_elections_ln",
                                      "Days_gov_approval_ln",
                                      "Draft_lenght_ln")]


colnames(COR_MATRIX) <- c("DV: Cabinet change", "Distance from coalition compromise (H1)","Saliency of a bill to other coalition parties (H2)", 
                           "Government's Legislative Plan (H3a)", "Ministerial phase change (H3b)",
                           "Involvement of LCG (H4)","Amendment", "Days to next parliamentary elections (N days logged)",
                           "Days in the cabinet (N logged)","Bill´s length (N words logged)")

rownames(COR_MATRIX) <-  c( "DV: Cabinet change", "Distance from coalition compromise (H1)","Saliency of a bill to other coalition parties (H2)", 
                            "Government's Legislative Plan (H3a)", "Ministerial phase change (H3b)",
                            "Involvement of LCG (H4)","Amendment", "Days to next parliamentary elections (N days logged)",
                            "Days in the cabinet (N logged)","Bill´s length (N words logged)")

COR_MATRIX <-cor(COR_MATRIX, use="complete.obs", method="spearman")
corrplot(COR_MATRIX, method="number", type="upper", tl.col = "black",  number.cex = 0.8, cl.cex = 0.9, tl.cex = 0.8)




#############################################################################################################
######################################### ANALYSIS    #######################################################
#############################################################################################################

################################################
# FIGURE 3 - Final model for the cabinet phase 
################################################

#  A new model after review:
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + SALIENCY_COAL_im_mean + Distance_COL_COMP_s + LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln + Legislative_plan +
                          as.factor(Government) +
                          as.factor(Proposer_cat) +
                          as.factor(PolicyArea),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)


library(coefplot)
coefplot::coefplot(Cabinet_phase_bi,
                   names=c( "Cabinet phase (bi-grams)"),  title = "DV: Draft change (0-100)",
                   legend.reverse = TRUE, intercept=FALSE,
                   xlab = "Standardized coefficients (Gelman method) + 90% c.i.", sort = c("magnitude"),
                   pointSize = 3, color = "black", zeroColor = "black",  innerCI = 1.645,   outerCI = 1.645,
                   coefficients=c(
                     "SALIENCY_COAL_im_mean",
                     "Distance_COL_COMP_s",
                     "Complexity", 
                     "Amended_law", 
                     "LRV_meeting",
                     "EU_transposition", 
                     "Days_to_elections_ln",
                     "Days_gov_approval_ln",
                     "Draft_lenght_ln",
                     "Legislative_plan",
                     "Draft_DIFF_BI"), 
                   newNames=c(Draft_DIFF_BI="Ministerial phase change (H3b)",
                              SALIENCY_COAL_im_mean = "Saliency of a bill to other coalition parties (H2)",
                              Distance_COL_COMP_s="Distance from coalition compromise (H1)",
                              LRV_meeting= "Involvement of LCG (H4)",
                              Amended_law="Amendment", 
                              Days_to_elections_ln="Days to next parliamentary elections (N days logged)",
                              Days_gov_approval_ln="Days in the cabinet (N logged)",
                              Draft_lenght_ln="Bill´s length (N words logged)",
                              Legislative_plan="Government's Legislative Plan (H3a)")) + theme_bw() +
  scale_shape_manual(values=c(17,19), guide = guide_legend(reverse = TRUE)) +
  theme(legend.position="bottom")  + 
  theme(plot.title = element_text(size=8,face="bold")) +
  theme(axis.text=element_text(size=8,face="bold"), axis.title=element_text(size=8,face="bold"))  +
  theme(legend.position="bottom") + scale_color_manual(values=c("black","black"), guide=FALSE)



# Comparing both phases - in the Online Appendix
Cabinet_phase <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + SALIENCY_COAL_im_mean + Distance_COL_COMP_s + LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln +
                          as.factor(Government) +
                          as.factor(Proposer_cat) +
                          as.factor(PolicyArea),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase)

multiplot(Ministerial_phase, Cabinet_phase,  
          names=c("Ministerial phase", "Cabinet phase"),  title = "DV: Draft change (0-100)",
          legend.reverse = TRUE, intercept=FALSE,
          xlab = "Standardized coefficients (Gelman method) + 90% c.i.", sort = c("magnitude"),
          pointSize = 3, color = "black", zeroColor = "black",  innerCI = 1.645,   outerCI = 1.645,
          coefficients=c( 
                         "SALIENCY_COAL_im_mean",
                         "Distance_COL_COMP_s",
                         "Complexity", 
                         "Amended_law", 
                         "LRV_meeting",
                         "EU_transposition", 
                         "Days_to_elections_ln",
                         "Days_gov_approval_ln",
                         "Draft_lenght_ln",
                         "Legislative_plan",
                         "Draft_DIFF_BI"), 
          newNames=c(Draft_DIFF_BI="Ministerial phase change (H3b)",
                     SALIENCY_COAL_im_mean = "Saliency of a bill to other coalition parties (H2)",
                     Distance_COL_COMP_s="Distance from coalition compromise (H1)",
                     LRV_meeting= "Involvement of LCG (H4)",
                     Amended_law="Amendment", 
                     Days_to_elections_ln="Days to next parliamentary elections (N days logged)",
                     Days_gov_approval_ln="Days in the cabinet (N logged)",
                     Draft_lenght_ln="Bill´s length (N words logged)",
                     Legislative_plan="Government's Legislative Plan (H3a)")) + theme_bw() +
  scale_shape_manual(values=c(17,19), guide = guide_legend(reverse = TRUE)) +
  theme(legend.position="bottom")  + 
  theme(plot.title = element_text(size=8,face="bold")) +
  theme(axis.text=element_text(size=8,face="bold"), axis.title=element_text(size=8,face="bold"))  +
  theme(legend.position="bottom") + scale_color_manual(values=c("black","black"), guide=FALSE)




# Testing control dummy variables - models with and without dummy variables
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + SALIENCY_COAL_im_mean + Distance_COL_COMP_s + LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln + Legislative_plan +
                          as.factor(Government) +
                          as.factor(Dimenze_Chapel_Hill_OTHER_CAT),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)

Cabinet_phase_bi_nodummy <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + SALIENCY_COAL_im_mean + Distance_COL_COMP_s + LRV_meeting +
                                  Amended_law + 
                                  Days_to_elections_ln +
                                  Days_gov_approval_ln +
                                  Draft_lenght_ln + Legislative_plan,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi_nodummy)

screenreg(list(Cabinet_phase_bi, Cabinet_phase_bi_nodummy),
          digits = 3,
          custom.model.names = c("FINAL", "no dummy"),
          single.row = TRUE) 


# Testing interaction effects:
# Figure 4 #
# Distance_COL_COMP_s*SALIENCY_COAL_im_mean
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + Saliency_CH_H + Distance_COL_COMP_s*SALIENCY_COAL_im_mean +  LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln
                          SOBOTKA_GOV + Dimenze_Chapel_Hill_1 + Dimenze_Chapel_Hill_2 + Dimenze_Chapel_Hill_3 + Dimenze_Chapel_Hill_4 + 
                          Dimenze_Chapel_Hill_5 + Dimenze_Chapel_Hill_6 + Dimenze_Chapel_Hill_7 + Dimenze_Chapel_Hill_8 + Dimenze_Chapel_Hill_9 +
                          Dimenze_Chapel_Hill_10,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)

interactions::interact_plot(Cabinet_phase_bi, pred = "Distance_COL_COMP_s", modx = "SALIENCY_COAL_im_mean",  modx.values = c(1, 3.3), modx.labels = c("min value", "max"),
                            interval= TRUE, legend.main = "Saliency (mean) for the coalition partners", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Distance from coalition compromise") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 


# Figure 5 # 
# Draft_DIFF_BI*SALIENCY_COAL_im_mean
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI*SALIENCY_COAL_im_mean + LRV_meeting + Saliency_CH_H + SALIENCY_COAL_im_mean  + Distance_COL_COMP_s +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln +
                          SOBOTKA_GOV + Dimenze_Chapel_Hill_1 + Dimenze_Chapel_Hill_2 + Dimenze_Chapel_Hill_3 + Dimenze_Chapel_Hill_4 + 
                          Dimenze_Chapel_Hill_5 + Dimenze_Chapel_Hill_6 + Dimenze_Chapel_Hill_7 + Dimenze_Chapel_Hill_8 + Dimenze_Chapel_Hill_9 +
                          Dimenze_Chapel_Hill_10,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)

interactions::interact_plot(Cabinet_phase_bi, pred = "Draft_DIFF_BI", modx = "SALIENCY_COAL_im_mean",  modx.values = c(1, 3.3), modx.labels = c("min value", "max"),
                            plot.points = FALSE,
                            interval= TRUE, legend.main = "Saliency (mean) for the coalition partners", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Change in ministerial phase") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 



##########################################
# Creating Table A - in the online appendix!
# Classical table format for regression result presentation

Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Distance_COL_COMP_s  + SALIENCY_COAL_im_mean + Legislative_plan + Draft_DIFF_BI  + LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln + 
                          as.factor(Government) +
                          as.factor(Proposer_cat) +
                          as.factor(PolicyArea),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)

library(texreg) 
htmlreg(list(Cabinet_phase_bi), # models
        file = "./FINAL_MODEL_TABLE_A_APPENDIX.html", # output file
        custom.model.names = c("Cabinet phase"), # model names
        reorder.coef = c(2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 
                         26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41,  1),
        custom.coef.names = c("(Constant)", "Distance from coalition compromise (H1)","Saliency of a bill to other coalition parties (H2)", 
                              "Government's Legislative Plan (H3a)", "Ministerial phase change (H3b)",
                              "Involvement of LCG (H4)","Amendment", "Days to next parliamentary elections (N days logged)",
                              "Days in the cabinet (N logged)","Bill´s length (N words logged)",
                              "VAR1", "VAR2", "VAR3","VAR4","VAR5","VAR6",
                              "VAR7","VAR8","VAR9","VAR10","VAR11","VAR12","VAR13","VAR14","VAR15",
                              "VAR16","VAR17","VAR18","VAR19","VAR20","VAR21","VAR22","VAR23","VAR24",
                              "VAR25","VAR26","VAR27","VAR28","VAR29","VAR30","VAR31"), # VAR - dummies
        custom.gof.names = c(" "," ", " ",  "Residual deviance","N"), # názvy proměnných
        digits = 3, # decimal places for coefficients and SEs
        single.row = FALSE, stars = c(0.01,  0.05, 0.1)) # display BETA and SEs under each other



# Regression diagnostic
##########################################

#VIF
car::vif(Cabinet_phase_bi)   #VIF is low

# Test of overdispersion
library(AER)
deviance(Cabinet_phase_bi)/Cabinet_phase_bi$df.residual
dispersiontest(Cabinet_phase_bi)


# Influencial plots
library(car)
influencePlot(Cabinet_phase_bi)

#Nonlinearity
crPlots(Cabinet_phase_bi)


# Normal Q-Q plot
res <- residuals(Cabinet_phase_bi, type="deviance")
plot(log(predict(Cabinet_phase_bi)), res)
abline(h=0, lty=2)
qqnorm(res)
qqline(res)

# general diagnostic plots
library(boot)
glm.diag.plots(Cabinet_phase_bi)



##### Testing other interactions:

# Saliency_CH_H*SALIENCY_COAL_im_mean 
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + Saliency_CH_H*SALIENCY_COAL_im_mean + Distance_COL_COMP_s + LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln +
                          SOBOTKA_GOV + Dimenze_Chapel_Hill_1 + Dimenze_Chapel_Hill_2 + Dimenze_Chapel_Hill_3 + Dimenze_Chapel_Hill_4 + 
                          Dimenze_Chapel_Hill_5 + Dimenze_Chapel_Hill_6 + Dimenze_Chapel_Hill_7 + Dimenze_Chapel_Hill_8 + Dimenze_Chapel_Hill_9 +
                          Dimenze_Chapel_Hill_10,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)


library(interactions)
interactions::interact_plot(Cabinet_phase_bi, pred = "Saliency_CH_H", modx = "SALIENCY_COAL_im_mean")
interactions::interact_plot(Cabinet_phase_bi, pred = "Saliency_CH_H", modx = "SALIENCY_COAL_im_mean", modx.values = c(1.3, 3.3), modx.labels = c("min value", "max"),
                            interval= TRUE, legend.main = "Saliency (mean) for the coalition partners", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Saliency of a draft for a drafting party") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 

summary(BILL_CHANGE_DATASET$SALIENCY_COAL_im_mean)
summary(BILL_CHANGE_DATASET$Saliency_CH_H)


interactions::interact_plot(Cabinet_phase_bi, pred = "SALIENCY_COAL_im_mean", modx = "Saliency_CH_H", modx.values = c(-1, 1), modx.labels = c("min value", "max"),
                            interval= TRUE, legend.main = "Saliency of a draft for a drafting party", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Saliency (mean) for the coalition partners") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 



# Saliency_CH_H*Distance_COL_COMP_s

Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + Saliency_CH_H*Distance_COL_COMP_s + SALIENCY_COAL_im_mean +  LRV_meeting +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln +
                          SOBOTKA_GOV + Dimenze_Chapel_Hill_1 + Dimenze_Chapel_Hill_2 + Dimenze_Chapel_Hill_3 + Dimenze_Chapel_Hill_4 + 
                          Dimenze_Chapel_Hill_5 + Dimenze_Chapel_Hill_6 + Dimenze_Chapel_Hill_7 + Dimenze_Chapel_Hill_8 + Dimenze_Chapel_Hill_9 +
                          Dimenze_Chapel_Hill_10,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)


library(interactions)
interactions::interact_plot(Cabinet_phase_bi, pred = "Distance_COL_COMP_s", modx = "Saliency_CH_H",  modx.values = c(-1.2, 1), modx.labels = c("min value", "max"),
                            interval= TRUE, legend.main = "Saliency of a draft for a drafting party", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Distance from coalition compromise") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 



# LRV_meeting*Distance_COL_COMP_s
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + Saliency_CH_H + SALIENCY_COAL_im_mean +  LRV_meeting*Distance_COL_COMP_s +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln +
                          SOBOTKA_GOV + Dimenze_Chapel_Hill_1 + Dimenze_Chapel_Hill_2 + Dimenze_Chapel_Hill_3 + Dimenze_Chapel_Hill_4 + 
                          Dimenze_Chapel_Hill_5 + Dimenze_Chapel_Hill_6 + Dimenze_Chapel_Hill_7 + Dimenze_Chapel_Hill_8 + Dimenze_Chapel_Hill_9 +
                          Dimenze_Chapel_Hill_10,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)

interactions::interact_plot(Cabinet_phase_bi, pred = "Distance_COL_COMP_s", modx = "LRV_meeting",  modx.values = c(-0.37, 1.76), modx.labels = c("min value", "max"),
                            interval= TRUE, legend.main = "Number of LRV meetings", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Distance from coalition compromise") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 

summary(BILL_CHANGE_DATASET$LRV_meeting)



# Saliency_CH_H*LRV_meeting
Cabinet_phase_bi <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Draft_DIFF_BI + Saliency_CH_H*LRV_meeting + SALIENCY_COAL_im_mean  + Distance_COL_COMP_s +
                          Amended_law + 
                          Days_to_elections_ln +
                          Days_gov_approval_ln +
                          Draft_lenght_ln +
                          SOBOTKA_GOV + Dimenze_Chapel_Hill_1 + Dimenze_Chapel_Hill_2 + Dimenze_Chapel_Hill_3 + Dimenze_Chapel_Hill_4 + 
                          Dimenze_Chapel_Hill_5 + Dimenze_Chapel_Hill_6 + Dimenze_Chapel_Hill_7 + Dimenze_Chapel_Hill_8 + Dimenze_Chapel_Hill_9 +
                          Dimenze_Chapel_Hill_10,  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_bi)

interactions::interact_plot(Cabinet_phase_bi, pred = "Saliency_CH_H", modx = "LRV_meeting",  modx.values = c(-0.37, 1.76), modx.labels = c("min value", "max"),
                            interval= TRUE, legend.main = "Number of LRV meetings", line.thickness = 0.8, colors = c("gray", "black")) + theme_gray() + 
  xlab("Saliency of a draft for a drafting party") + ylab("Predicted change in the cabinet phase") + 
  theme(plot.title = element_text(size=16,face="bold")) +
  theme(axis.text=element_text(size=13), axis.title=element_text(size=13,face="bold")) + 
  theme(legend.position="top") 





# Predicted values based on values of X controling for Z
########################################################

# Cabinet phase
dat <- ggpredict(model = Cabinet_phase,
                 terms = "Distance_COL_COMP_s ",
                 ci.lvl = 0.95)
plot(dat)
dat <- ggpredict(model = Cabinet_phase,
                 terms = "SALIENCY_COAL_im_mean",
                 ci.lvl = 0.95)
plot(dat)



dat <- ggpredict(model = Cabinet_phase,
                 terms = "Legislative_plan",
                 ci.lvl = 0.95)
plot(dat)

dat <- ggpredict(model = Cabinet_phase,
                 terms = "Draft_DIFF_BI  ",
                 ci.lvl = 0.95)
plot(dat)

dat <- ggpredict(model = Cabinet_phase,
                 terms = "LRV_meeting",
                 ci.lvl = 0.95)
plot(dat)



dat <- ggpredict(model = Cabinet_phase,
                 terms = "Amended_law",
                 ci.lvl = 0.95)
plot(dat)



dat <- ggpredict(model = Cabinet_phase,
                 terms = "Days_to_elections_ln",
                 ci.lvl = 0.95)
plot(dat)

dat <- ggpredict(model = Cabinet_phase,
                 terms = "Days_gov_approval_ln",
                 ci.lvl = 0.95)
plot(dat)

dat <- ggpredict(model = Cabinet_phase,
                 terms = "Draft_lenght_ln",
                 ci.lvl = 0.95)
plot(dat)




# Descriptive means plots with confidence intervals
######################################################################################################################

#Run function before! - computes standard errors
stderr <- function(x, na.rm=TRUE) {
  if (na.rm) x <- na.omit(x)
  sqrt(var(x)/length(x))
}



#Czech abbreviations of ministries
BILL_CHANGE_DATASET$Proposer_cat <- factor(BILL_CHANGE_DATASET$Proposer_cat,
                            levels = c(2,3,4,5,6,7,8,9,10,11,12,14,15),
                            labels = c("MD", "MF", "MK", "MMR", "MO", "MPO", "MPSV", "MSP", "MV", "MZD", "MZE", "MŠMT", "MŽP"))
#English labels of ministries

BILL_CHANGE_DATASET$Proposer_cat <- factor(BILL_CHANGE_DATASET$Proposer_cat,
                                           levels = c(2,3,4,5,6,7,8,9,10,11,12,14,15),
                                           labels = c("Transportation", "Finance", "Culture", "Regional Dev.", "Defence",
                                                      "Industry", "Labour Social Aff.", "Justice", "Interior", "Health",
                                                      "Agriculture", "Education", "Environment"))


# Draft change - Ministerial phase by ministry

BILL_CHANGE_DATASET %>%
  group_by(Proposer_cat) %>%
  summarise(EFF = mean(Draft_DIFF, na.rm = TRUE),
            SE = stderr(Draft_DIFF)) %>%
  ggplot(aes(x = reorder(Proposer_cat, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Ministry",
       y = "Draft change - Ministerial phase") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))

# Draft change - Cabinet phase by ministry
BILL_CHANGE_DATASET %>%
  group_by(Proposer_cat) %>%
  summarise(EFF = mean(DRAFT_DIFF_GOV_PAR, na.rm = TRUE),
            SE = stderr(DRAFT_DIFF_GOV_PAR)) %>%
  ggplot(aes(x = reorder(Proposer_cat, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Ministry",
       y = "Draft change - Cabinet phase") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))

# Draft change - Executive phase (total) by ministry
BILL_CHANGE_DATASET %>%
  group_by(Proposer_cat) %>%
  summarise(EFF = mean(DRAFT_DIFF_IMC_PAR, na.rm = TRUE),
            SE = stderr(DRAFT_DIFF_IMC_PAR)) %>%
  ggplot(aes(x = reorder(Proposer_cat, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Ministry",
       y = "Draft change - Executive phase (total)") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))






# Distance from coalition compromise -  by ministries (drafters)
BILL_CHANGE_DATASET %>%
  group_by(Proposer_cat) %>%
  summarise(EFF = mean(Distance_COL_COMP_s, na.rm = TRUE),
            SE = stderr(Distance_COL_COMP_s)) %>%
  ggplot(aes(x = reorder(Proposer_cat, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Ministry",
       y = "Distance from coalition compromise") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))



# Draft change -  by Chapel hill dimensions
BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT <- factor(BILL_CHANGE_DATASET$Dimenze_Chapel_Hill_OTHER_CAT,
                                                            levels = c(2,3,4,5,6,7,8,9,10),
                                                            labels = c("civlib_laworder", "deregulation", "environment", "other",
                                                                       "redistribution", "regions", "sociallifestyle", "spendvtax", "urban_rural"))


BILL_CHANGE_DATASET %>%
  group_by(Dimenze_Chapel_Hill_OTHER_CAT) %>%
  summarise(EFF = mean(Draft_DIFF, na.rm = TRUE),
            SE = stderr(Draft_DIFF)) %>%
  ggplot(aes(x = reorder(Dimenze_Chapel_Hill_OTHER_CAT, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Chapel Hill dimensions",
       y = "Draft change - Ministerial phase") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))

BILL_CHANGE_DATASET %>%
  group_by(Dimenze_Chapel_Hill_OTHER_CAT) %>%
  summarise(EFF = mean(DRAFT_DIFF_GOV_PAR, na.rm = TRUE),
            SE = stderr(DRAFT_DIFF_GOV_PAR)) %>%
  ggplot(aes(x = reorder(Dimenze_Chapel_Hill_OTHER_CAT, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Chapel Hill dimensions",
       y = "Draft change - Cabinet phase") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))

BILL_CHANGE_DATASET %>%
  group_by(Dimenze_Chapel_Hill_OTHER_CAT) %>%
  summarise(EFF = mean(DRAFT_DIFF_IMC_PAR, na.rm = TRUE),
            SE = stderr(DRAFT_DIFF_IMC_PAR)) %>%
  ggplot(aes(x = reorder(Dimenze_Chapel_Hill_OTHER_CAT, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Chapel Hill dimensions",
       y = "Draft change - Executive phase (total)") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))
 

# Distance from coalition compromise -  by Chapel hill dimensions
BILL_CHANGE_DATASET %>%
  group_by(Dimenze_Chapel_Hill_OTHER_CAT) %>%
  summarise(EFF = mean(Distance_COL_COMP_s, na.rm = TRUE),
            SE = stderr(Distance_COL_COMP_s)) %>%
  ggplot(aes(x = reorder(Dimenze_Chapel_Hill_OTHER_CAT, -EFF),
             y = EFF)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = EFF - 1.96*SE,
                    ymax = EFF + 1.96*SE),
                width = 0, size = 1.25) +
  labs(x = "Chapel Hill dimensions",
       y = "Distance from coalition compromise") +
  theme(plot.title = element_text(size=14,face="bold")) +              
  theme(axis.text=element_text(size=10, angle=45, hjust = 1), axis.title=element_text(size=13,face="bold"))




# Cabinet phase

Cabinet_phase_UNI <- glm(DRAFT_DIFF_GOV_PAR_UNI ~ Main_ruling_party_DB +
                       Comments_coal_partner_ln +
                       Comments_same_party_ln + 
                       COM_IMP_DB +
                       Saliency_CH_H +
                       Distance_COL_COMP_s +
                       Complexity + 
                       Amended_law + 
                       LRV_meeting +
                       EU_transposition + 
                       Days_to_elections_ln +
                       Days_gov_approval_ln +
                       Draft_lenght_ln +
                       Legislative_plan + 
                       as.factor(Government) + as.factor(Proposer_cat) +
                       as.factor(Dimenze_Chapel_Hill_OTHER_CAT) + as.factor(PolicyArea),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_UNI)


Cabinet_phase_BI <- glm(DRAFT_DIFF_GOV_PAR_BI ~ Main_ruling_party_DB +
                       Comments_coal_partner_ln +
                       Comments_same_party_ln + 
                       COM_IMP_DB +
                       Saliency_CH_H +
                       Distance_COL_COMP_s +
                       Complexity + 
                       Amended_law + 
                       LRV_meeting +
                       EU_transposition + 
                       Days_to_elections_ln +
                       Days_gov_approval_ln +
                       Draft_lenght_ln +
                       Legislative_plan + 
                       as.factor(Government) + as.factor(Proposer_cat) +
                       as.factor(Dimenze_Chapel_Hill_OTHER_CAT) + as.factor(PolicyArea),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase_BI)


Cabinet_phase <- glm(DRAFT_DIFF_GOV_PAR ~ Main_ruling_party_DB +
                       Comments_coal_partner_ln +
                       Comments_same_party_ln + 
                       COM_IMP_DB +
                       Saliency_CH_H +
                       Distance_COL_COMP_s +
                       Complexity + 
                       Amended_law + 
                       LRV_meeting +
                       EU_transposition + 
                       Days_to_elections_ln +
                       Days_gov_approval_ln +
                       Draft_lenght_ln +
                       Legislative_plan + 
                       as.factor(Government) + as.factor(Proposer_cat) +
                       as.factor(Dimenze_Chapel_Hill_OTHER_CAT) + as.factor(PolicyArea),  family=poisson(link=log), data=BILL_CHANGE_DATASET)
summary(Cabinet_phase)



screenreg(list(Cabinet_phase_UNI, Cabinet_phase_BI, Cabinet_phase),
          digits = 3,
          custom.model.names = c("Cabinet phase - unigrams", "Cabinet phase - bigrams" , "Cabinet phase - starý"),
          single.row = TRUE) # See change when setting this to FALSE


